{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# Load all the required packages\n",
    "using ODEInterface\n",
    "using ForwardDiff\n",
    "using Gadfly\n",
    "using Colors\n",
    "@ODEInterface.import_huge\n",
    "loadODESolvers();\n",
    "\n",
    "# Define the right-hand function for automatic differentiation\n",
    "function vdpolAD(x)\n",
    "    return [x[2],((1-x[1]^2)*x[2]-x[1])*1e6]\n",
    "end\n",
    "\n",
    "# Define the system for the solver\n",
    "function vdpol(t,x,dx)\n",
    "    dx[:] = vdpolAD(x);\n",
    "    return nothing\n",
    "end\n",
    "\n",
    "# Define the Jacobian function using AD\n",
    "function getJacobian(t,x,J)\n",
    "    J[:,:] = ForwardDiff.jacobian(vdpolAD,x);\n",
    "    return nothing\n",
    "end\n",
    "\n",
    "# Flag to check whether plot is to be generated and saved or not\n",
    "# Also checks if all solvers are successful\n",
    "printFlag = true;\n",
    "\n",
    "# Initial conditions\n",
    "t0 = 0.0; T = [1.0:11.0;]; x0 = [2.0,0.0];\n",
    "\n",
    "# Get \"reference solution\"\n",
    "Tol = 1e-14;\n",
    "# for Tol < 1e-14 we get the error \"TOLERANCES ARE TOO SMALL\"\n",
    "opt = OptionsODE(OPT_EPS=>1.11e-16,OPT_RTOL=>Tol, OPT_ATOL=>Tol,\n",
    "OPT_RHS_CALLMODE => RHS_CALL_INSITU,\n",
    "OPT_JACOBIMATRIX => getJacobian);\n",
    "\n",
    "(t,x,retcode,stats) = odecall(seulex,vdpol,[t0, T[end]], x0, opt);\n",
    "t = [t0;t];\n",
    "x = [x0';x];"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "p = plot(x=t,y=x[:,1],Geom.path,\n",
    "Theme(line_width=2pt,default_color=colorant\"black\",\n",
    "panel_stroke=colorant\"black\",key_position=:top,\n",
    "key_max_columns = 1,major_label_font_size=24pt,minor_label_font_size=22pt,\n",
    "key_title_font_size=22pt,key_label_font_size=20pt),\n",
    "Coord.cartesian(ymin=-2.5,ymax=2.5,xmin=0,xmax=11),\n",
    "Guide.xlabel(\"time (s)\"),Guide.ylabel(\"Position (x)\"),\n",
    "Guide.xticks(ticks=[0:11;]),Guide.yticks(ticks=[-2:0.5:2;]));\n",
    "\n",
    "draw(PNG(\"../../ImagesAndPDFs/Plots/vdpolPlot.png\",30cm,20cm),p);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Julia 0.4.2",
   "language": "julia",
   "name": "julia-0.4"
  },
  "language_info": {
   "file_extension": ".jl",
   "mimetype": "application/julia",
   "name": "julia",
   "version": "0.4.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}
